Load Data
dataset <- read.delim("raw_data/FigureS5N.txt", stringsAsFactors = FALSE)
dataset$genotype <- gsub(" ", "", dataset$genotype )
dataset$genotype <- factor(dataset$genotype)
dataset$Experiment <- factor(rep(paste0("exp", 1:(length(dataset$genotype)/length(levels(dataset$genotype)))), each=length(unique(dataset$genotype))))
dataset$UID <- factor(paste(dataset$Experiment, dataset$genotype))
# wide format
kable(dataset, row.names = F)
| WT |
1097 |
575 |
397 |
106 |
exp1 |
exp1 WT |
| YFP-ALC1 |
1425 |
1018 |
740 |
100 |
exp1 |
exp1 YFP-ALC1 |
| WT |
1768 |
942 |
462 |
95 |
exp2 |
exp2 WT |
| YFP-ALC1 |
1271 |
931 |
528 |
150 |
exp2 |
exp2 YFP-ALC1 |
| WT |
1252 |
741 |
456 |
95 |
exp3 |
exp3 WT |
| YFP-ALC1 |
1354 |
935 |
624 |
48 |
exp3 |
exp3 YFP-ALC1 |
library(reshape2)
# reshape to long format
dataset <- melt(dataset, variable.name = "Treatment", value.name = "Counts")
dataset$genotype <- relevel(dataset$genotype, ref = "WT")
dataset$UID <- relevel(dataset$UID, ref = "exp1 WT")
dataset$H2O2 <- gsub("NT","1",dataset$Treatment)
dataset$H2O2 <- gsub("H2O2_|uM","",dataset$H2O2)
dataset$H2O2 <- log10(as.integer(dataset$H2O2))
dataset$Offset <- NA
for(uid in levels(dataset$UID)){
dataset$Offset[dataset$UID == uid] <- mean(dataset$Counts[dataset$UID == uid])
}
dataset$NormCounts <- dataset$Counts / dataset$Offset
dataset$Offset2 <- NA
for(gid in levels(dataset$genotype)){
dataset$Offset2[dataset$genotype == gid] <- mean(dataset$NormCounts[dataset$genotype == gid & dataset$H2O2 == 0])
}
dataset$NormCounts2 <- dataset$NormCounts / dataset$Offset2
# long format
kable(dataset, row.names = F)
| WT |
exp1 |
exp1 WT |
NT |
1097 |
0.000000 |
543.75 |
2.0174713 |
2.050234 |
0.9840200 |
| YFP-ALC1 |
exp1 |
exp1 YFP-ALC1 |
NT |
1425 |
0.000000 |
820.75 |
1.7362169 |
1.776869 |
0.9771216 |
| WT |
exp2 |
exp2 WT |
NT |
1768 |
0.000000 |
816.75 |
2.1646771 |
2.050234 |
1.0558196 |
| YFP-ALC1 |
exp2 |
exp2 YFP-ALC1 |
NT |
1271 |
0.000000 |
720.00 |
1.7652778 |
1.776869 |
0.9934767 |
| WT |
exp3 |
exp3 WT |
NT |
1252 |
0.000000 |
636.00 |
1.9685535 |
2.050234 |
0.9601604 |
| YFP-ALC1 |
exp3 |
exp3 YFP-ALC1 |
NT |
1354 |
0.000000 |
740.25 |
1.8291118 |
1.776869 |
1.0294017 |
| WT |
exp1 |
exp1 WT |
H2O2_5uM |
575 |
0.698970 |
543.75 |
1.0574713 |
2.050234 |
0.5157808 |
| YFP-ALC1 |
exp1 |
exp1 YFP-ALC1 |
H2O2_5uM |
1018 |
0.698970 |
820.75 |
1.2403290 |
1.776869 |
0.6980419 |
| WT |
exp2 |
exp2 WT |
H2O2_5uM |
942 |
0.698970 |
816.75 |
1.1533517 |
2.050234 |
0.5625464 |
| YFP-ALC1 |
exp2 |
exp2 YFP-ALC1 |
H2O2_5uM |
931 |
0.698970 |
720.00 |
1.2930556 |
1.776869 |
0.7277158 |
| WT |
exp3 |
exp3 WT |
H2O2_5uM |
741 |
0.698970 |
636.00 |
1.1650943 |
2.050234 |
0.5682739 |
| YFP-ALC1 |
exp3 |
exp3 YFP-ALC1 |
H2O2_5uM |
935 |
0.698970 |
740.25 |
1.2630868 |
1.776869 |
0.7108498 |
| WT |
exp1 |
exp1 WT |
H2O2_30uM |
397 |
1.477121 |
543.75 |
0.7301149 |
2.050234 |
0.3561130 |
| YFP-ALC1 |
exp1 |
exp1 YFP-ALC1 |
H2O2_30uM |
740 |
1.477121 |
820.75 |
0.9016144 |
1.776869 |
0.5074175 |
| WT |
exp2 |
exp2 WT |
H2O2_30uM |
462 |
1.477121 |
816.75 |
0.5656566 |
2.050234 |
0.2758985 |
| YFP-ALC1 |
exp2 |
exp2 YFP-ALC1 |
H2O2_30uM |
528 |
1.477121 |
720.00 |
0.7333333 |
1.776869 |
0.4127110 |
| WT |
exp3 |
exp3 WT |
H2O2_30uM |
456 |
1.477121 |
636.00 |
0.7169811 |
2.050234 |
0.3497070 |
| YFP-ALC1 |
exp3 |
exp3 YFP-ALC1 |
H2O2_30uM |
624 |
1.477121 |
740.25 |
0.8429585 |
1.776869 |
0.4744067 |
| WT |
exp1 |
exp1 WT |
H2O2_50uM |
106 |
1.698970 |
543.75 |
0.1949425 |
2.050234 |
0.0950831 |
| YFP-ALC1 |
exp1 |
exp1 YFP-ALC1 |
H2O2_50uM |
100 |
1.698970 |
820.75 |
0.1218398 |
1.776869 |
0.0685699 |
| WT |
exp2 |
exp2 WT |
H2O2_50uM |
95 |
1.698970 |
816.75 |
0.1163147 |
2.050234 |
0.0567324 |
| YFP-ALC1 |
exp2 |
exp2 YFP-ALC1 |
H2O2_50uM |
150 |
1.698970 |
720.00 |
0.2083333 |
1.776869 |
0.1172474 |
| WT |
exp3 |
exp3 WT |
H2O2_50uM |
95 |
1.698970 |
636.00 |
0.1493711 |
2.050234 |
0.0728556 |
| YFP-ALC1 |
exp3 |
exp3 YFP-ALC1 |
H2O2_50uM |
48 |
1.698970 |
740.25 |
0.0648430 |
1.776869 |
0.0364928 |
Plot Data
library(ggplot2)
# raw data
ggplot(dataset, aes(x=H2O2, y=Counts)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE, aes(colour=genotype)) +
geom_point(aes(colour=genotype, shape=Experiment), size=2) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)") +
scale_shape_manual(values=15:20) +
scale_color_manual(values=c("#000000","#808000"))

# NormCounts Linear
ggplot(dataset, aes(x=H2O2, y=NormCounts, color=genotype)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=genotype), size=2) +
geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)") +
scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Linear
ggplot(dataset, aes(x=H2O2, y=NormCounts2, color=genotype)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=genotype), size=2) +
geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)") +
scale_color_manual(values=c("#000000","#808000"))

# NormCounts Quadratic
ggplot(dataset, aes(x=H2O2, y=NormCounts, color=genotype)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=genotype), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)")+
scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Quadratic
ggplot(dataset, aes(x=H2O2, y=NormCounts2, color=genotype)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=genotype), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)") +
scale_color_manual(values=c("#000000","#808000"))

# NormCounts Cubic
ggplot(dataset, aes(x=H2O2, y=NormCounts, color=genotype)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=genotype), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)")+
scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Cubic
ggplot(dataset, aes(x=H2O2, y=NormCounts2, color=genotype)) +
theme_bw() +
theme(panel.grid=element_blank(), text = element_text(size=14)) +
geom_point(aes(colour=genotype), size=2) +
geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)") +
scale_color_manual(values=c("#000000","#808000"))

library(Cairo)
cairo_pdf("FigureS5N.pdf", width = 5, height = 4, family = "Arial")
ggplot(dataset, aes(x=H2O2, y=NormCounts2)) +
theme_bw() +
theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(),
axis.line = element_line(colour = "black"), text = element_text(size=14),
panel.border = element_blank(), panel.background = element_blank()) +
geom_point(aes(colour = genotype)) +
geom_smooth(method=lm, formula = y ~ poly(x,3), se=TRUE, aes(colour = genotype), fill='#CCCCCC') +
#facet_grid(. ~ genotype) +
xlab(label = "H2O2 (log10 µM)") +
ylab(label = "Normalized Counts") +
scale_color_manual(values=c("#000000","#808000"))
dev.off()
## quartz_off_screen
## 2
Models
library(MASS)
library(DHARMa)
library(lme4)
library(lmerTest)
library(bbmle)
Linear formula
fit1 <- lm(Counts ~ Experiment + H2O2*genotype, data = dataset)
print(summary(fit1))
##
## Call:
## lm(formula = Counts ~ Experiment + H2O2 * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -245.04 -132.06 -19.98 95.30 382.01
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1299.866 110.649 11.748 7.10e-10 ***
## Experimentexp2 86.125 94.856 0.908 0.376
## Experimentexp3 5.875 94.856 0.062 0.951
## H2O2 -686.474 81.573 -8.416 1.18e-07 ***
## genotypeYFP-ALC1 60.660 135.971 0.446 0.661
## H2O2:genotypeYFP-ALC1 35.276 115.361 0.306 0.763
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 189.7 on 18 degrees of freedom
## Multiple R-squared: 0.8839, Adjusted R-squared: 0.8517
## F-statistic: 27.42 on 5 and 18 DF, p-value: 7.994e-08
cat("AIC: ", AIC(fit1))
## AIC: 326.9888
simres <- simulateResiduals(fittedModel = fit1)
plot(simres)

fit2 <- lm(NormCounts ~ H2O2*genotype, data = dataset)
print(summary(fit2))
##
## Call:
## lm(formula = NormCounts ~ H2O2 * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.30870 -0.11268 -0.01016 0.10500 0.33774
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.98687 0.09414 21.105 3.86e-15 ***
## H2O2 -1.01869 0.07987 -12.754 4.60e-11 ***
## genotypeYFP-ALC1 -0.15575 0.13314 -1.170 0.256
## H2O2:genotypeYFP-ALC1 0.16077 0.11295 1.423 0.170
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1858 on 20 degrees of freedom
## Multiple R-squared: 0.9329, Adjusted R-squared: 0.9228
## F-statistic: 92.68 on 3 and 20 DF, p-value: 6.638e-12
cat("AIC: ", AIC(fit2))
## AIC: -7.066552
simres <- simulateResiduals(fittedModel = fit2)
plot(simres)

fit3 <- lm(NormCounts2 ~ H2O2*genotype, data = dataset)
print(summary(fit3))
##
## Call:
## lm(formula = NormCounts2 ~ H2O2 * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.173733 -0.054959 -0.005033 0.052231 0.190077
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.96909 0.05038 19.234 2.27e-14 ***
## H2O2 -0.49686 0.04275 -11.624 2.38e-10 ***
## genotypeYFP-ALC1 0.06144 0.07125 0.862 0.399
## H2O2:genotypeYFP-ALC1 0.01404 0.06045 0.232 0.819
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.09941 on 20 degrees of freedom
## Multiple R-squared: 0.9301, Adjusted R-squared: 0.9196
## F-statistic: 88.7 on 3 and 20 DF, p-value: 9.977e-12
cat("AIC: ", AIC(fit3))
## AIC: -37.07284
simres <- simulateResiduals(fittedModel = fit3)
plot(simres)

fit4 <- lmer(Counts ~ H2O2*genotype + (1|UID), data = dataset)
print(summary(fit4))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ H2O2 * genotype + (1 | UID)
## Data: dataset
##
## REML criterion at convergence: 273.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.35434 -0.61511 -0.05378 0.58692 2.21547
##
## Random effects:
## Groups Name Variance Std.Dev.
## UID (Intercept) 3127 55.92
## Residual 31742 178.16
## Number of obs: 24, groups: UID, 6
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 1330.53 95.89 15.97 13.875 2.51e-10 ***
## H2O2 -686.47 76.61 16.00 -8.961 1.24e-07 ***
## genotypeYFP-ALC1 60.66 135.61 15.97 0.447 0.661
## H2O2:genotypeYFP-ALC1 35.28 108.34 16.00 0.326 0.749
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) H2O2 gYFP-A
## H2O2 -0.774
## gntYFP-ALC1 -0.707 0.547
## H2O2:YFP-AL 0.547 -0.707 -0.774
cat("AIC: ", AIC(fit4))
## AIC: 285.7402
simres <- simulateResiduals(fittedModel = fit4)
plot(simres)

Quadratic formula
fit5 <- lm(Counts ~ Experiment + poly(H2O2,2)*genotype, data = dataset)
print(summary(fit5))
##
## Call:
## lm(formula = Counts ~ Experiment + poly(H2O2, 2) * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -226.0 -113.0 -48.0 102.0 358.9
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 634.833 77.540 8.187 4.1e-07 ***
## Experimentexp2 86.125 94.966 0.907 0.378
## Experimentexp3 5.875 94.966 0.062 0.951
## poly(H2O2, 2)1 -2257.814 268.605 -8.406 2.9e-07 ***
## poly(H2O2, 2)2 121.603 268.605 0.453 0.657
## genotypeYFP-ALC1 94.833 77.540 1.223 0.239
## poly(H2O2, 2)1:genotypeYFP-ALC1 116.021 379.865 0.305 0.764
## poly(H2O2, 2)2:genotypeYFP-ALC1 -477.253 379.865 -1.256 0.227
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 189.9 on 16 degrees of freedom
## Multiple R-squared: 0.8966, Adjusted R-squared: 0.8513
## F-statistic: 19.82 on 7 and 16 DF, p-value: 9.13e-07
cat("AIC: ", AIC(fit5))
## AIC: 328.218
simres <- simulateResiduals(fittedModel = fit5)
plot(simres)

fit6 <- lm(NormCounts ~ poly(H2O2,2)*genotype, data = dataset)
print(summary(fit6))
##
## Call:
## lm(formula = NormCounts ~ poly(H2O2, 2) * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.22690 -0.11902 -0.05308 0.10590 0.31224
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.000e+00 5.115e-02 19.552 1.42e-13 ***
## poly(H2O2, 2)1 -3.350e+00 2.506e-01 -13.372 8.67e-11 ***
## poly(H2O2, 2)2 1.813e-01 2.506e-01 0.724 0.4786
## genotypeYFP-ALC1 2.783e-16 7.233e-02 0.000 1.0000
## poly(H2O2, 2)1:genotypeYFP-ALC1 5.288e-01 3.543e-01 1.492 0.1530
## poly(H2O2, 2)2:genotypeYFP-ALC1 -6.474e-01 3.543e-01 -1.827 0.0843 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1772 on 18 degrees of freedom
## Multiple R-squared: 0.9451, Adjusted R-squared: 0.9298
## F-statistic: 61.92 on 5 and 18 DF, p-value: 1.043e-10
cat("AIC: ", AIC(fit6))
## AIC: -7.865543
simres <- simulateResiduals(fittedModel = fit6)
plot(simres)

fit7 <- lm(NormCounts2 ~ poly(H2O2,2)*genotype, data = dataset)
print(summary(fit7))
##
## Call:
## lm(formula = NormCounts2 ~ poly(H2O2, 2) * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.12770 -0.06473 -0.02589 0.05405 0.17573
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.48775 0.02716 17.958 6.13e-13 ***
## poly(H2O2, 2)1 -1.63419 0.13306 -12.282 3.47e-10 ***
## poly(H2O2, 2)2 0.08843 0.13306 0.665 0.5147
## genotypeYFP-ALC1 0.07504 0.03841 1.954 0.0665 .
## poly(H2O2, 2)1:genotypeYFP-ALC1 0.04617 0.18818 0.245 0.8089
## poly(H2O2, 2)2:genotypeYFP-ALC1 -0.35076 0.18818 -1.864 0.0787 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.09409 on 18 degrees of freedom
## Multiple R-squared: 0.9436, Adjusted R-squared: 0.928
## F-statistic: 60.28 on 5 and 18 DF, p-value: 1.308e-10
cat("AIC: ", AIC(fit7))
## AIC: -38.24457
simres <- simulateResiduals(fittedModel = fit7)
plot(simres)

fit8 <- lmer(Counts ~ poly(H2O2,2)*genotype + (1|UID), data = dataset)
print(summary(fit8))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(H2O2, 2) * genotype + (1 | UID)
## Data: dataset
##
## REML criterion at convergence: 241
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.2496 -0.5899 -0.2551 0.5390 2.0928
##
## Random effects:
## Groups Name Variance Std.Dev.
## UID (Intercept) 3255 57.05
## Residual 31231 176.72
## Number of obs: 24, groups: UID, 6
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 665.50 60.73 4.00 10.959 0.000394
## poly(H2O2, 2)1 -2257.81 249.92 14.00 -9.034 3.24e-07
## poly(H2O2, 2)2 121.60 249.92 14.00 0.487 0.634098
## genotypeYFP-ALC1 94.83 85.88 4.00 1.104 0.331429
## poly(H2O2, 2)1:genotypeYFP-ALC1 116.02 353.45 14.00 0.328 0.747575
## poly(H2O2, 2)2:genotypeYFP-ALC1 -477.25 353.45 14.00 -1.350 0.198350
##
## (Intercept) ***
## poly(H2O2, 2)1 ***
## poly(H2O2, 2)2
## genotypeYFP-ALC1
## poly(H2O2, 2)1:genotypeYFP-ALC1
## poly(H2O2, 2)2:genotypeYFP-ALC1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) pl(H2O2,2)1 pl(H2O2,2)2 gYFP-A p(H2O2,2)1:
## pl(H2O2,2)1 0.000
## pl(H2O2,2)2 0.000 0.000
## gntYFP-ALC1 -0.707 0.000 0.000
## p(H2O2,2)1: 0.000 -0.707 0.000 0.000
## p(H2O2,2)2: 0.000 0.000 -0.707 0.000 0.000
cat("AIC: ", AIC(fit8))
## AIC: 256.9578
simres <- simulateResiduals(fittedModel = fit8)
plot(simres)

Cubic formula
fit9 <- lm(Counts ~ Experiment + poly(H2O2,3)*genotype, data = dataset)
print(summary(fit9))
##
## Call:
## lm(formula = Counts ~ Experiment + poly(H2O2, 3) * genotype,
## data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -244.67 -65.79 5.79 39.11 340.21
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 634.833 62.149 10.215 7.18e-08 ***
## Experimentexp2 86.125 76.117 1.131 0.2769
## Experimentexp3 5.875 76.117 0.077 0.9396
## poly(H2O2, 3)1 -2257.814 215.292 -10.487 5.17e-08 ***
## poly(H2O2, 3)2 121.603 215.292 0.565 0.5811
## poly(H2O2, 3)3 -410.102 215.292 -1.905 0.0775 .
## genotypeYFP-ALC1 94.833 62.149 1.526 0.1493
## poly(H2O2, 3)1:genotypeYFP-ALC1 116.021 304.468 0.381 0.7089
## poly(H2O2, 3)2:genotypeYFP-ALC1 -477.253 304.468 -1.567 0.1393
## poly(H2O2, 3)3:genotypeYFP-ALC1 -170.664 304.468 -0.561 0.5840
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 152.2 on 14 degrees of freedom
## Multiple R-squared: 0.9419, Adjusted R-squared: 0.9045
## F-statistic: 25.2 on 9 and 14 DF, p-value: 3.794e-07
cat("AIC: ", AIC(fit9))
## AIC: 318.3933
simres <- simulateResiduals(fittedModel = fit9)
plot(simres)

fit10 <- lm(NormCounts ~ poly(H2O2,3)*genotype, data = dataset)
print(summary(fit10))
##
## Call:
## lm(formula = NormCounts ~ poly(H2O2, 3) * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.105261 -0.038084 -0.003288 0.042566 0.114443
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.000e+00 2.023e-02 49.429 < 2e-16 ***
## poly(H2O2, 3)1 -3.350e+00 9.911e-02 -33.805 2.61e-16 ***
## poly(H2O2, 3)2 1.813e-01 9.911e-02 1.829 0.086062 .
## poly(H2O2, 3)3 -6.343e-01 9.911e-02 -6.400 8.79e-06 ***
## genotypeYFP-ALC1 -1.694e-16 2.861e-02 0.000 1.000000
## poly(H2O2, 3)1:genotypeYFP-ALC1 5.288e-01 1.402e-01 3.773 0.001667 **
## poly(H2O2, 3)2:genotypeYFP-ALC1 -6.474e-01 1.402e-01 -4.619 0.000284 ***
## poly(H2O2, 3)3:genotypeYFP-ALC1 -1.210e-01 1.402e-01 -0.863 0.400712
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.07008 on 16 degrees of freedom
## Multiple R-squared: 0.9924, Adjusted R-squared: 0.989
## F-statistic: 296.8 on 7 and 16 DF, p-value: 1.015e-15
cat("AIC: ", AIC(fit10))
## AIC: -51.21038
simres <- simulateResiduals(fittedModel = fit10)
plot(simres)

fit11 <- lm(NormCounts2 ~ poly(H2O2,3)*genotype, data = dataset)
print(summary(fit11))
##
## Call:
## lm(formula = NormCounts2 ~ poly(H2O2, 3) * genotype, data = dataset)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.052134 -0.019338 -0.001694 0.020761 0.055820
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.48775 0.01049 46.489 < 2e-16 ***
## poly(H2O2, 3)1 -1.63419 0.05140 -31.794 6.87e-16 ***
## poly(H2O2, 3)2 0.08843 0.05140 1.720 0.104632
## poly(H2O2, 3)3 -0.30939 0.05140 -6.019 1.79e-05 ***
## genotypeYFP-ALC1 0.07504 0.01484 5.057 0.000117 ***
## poly(H2O2, 3)1:genotypeYFP-ALC1 0.04617 0.07269 0.635 0.534263
## poly(H2O2, 3)2:genotypeYFP-ALC1 -0.35076 0.07269 -4.825 0.000186 ***
## poly(H2O2, 3)3:genotypeYFP-ALC1 -0.11570 0.07269 -1.592 0.131009
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.03634 on 16 degrees of freedom
## Multiple R-squared: 0.9925, Adjusted R-squared: 0.9893
## F-statistic: 303.5 on 7 and 16 DF, p-value: 8.507e-16
cat("AIC: ", AIC(fit11))
## AIC: -82.72848
simres <- simulateResiduals(fittedModel = fit11)
plot(simres)

fit12 <- lmer(Counts ~ poly(H2O2,3)*genotype + (1|UID), data = dataset)
print(summary(fit12))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(H2O2, 3) * genotype + (1 | UID)
## Data: dataset
##
## REML criterion at convergence: 204.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.57978 -0.46992 0.05654 0.30014 2.39500
##
## Random effects:
## Groups Name Variance Std.Dev.
## UID (Intercept) 7219 84.96
## Residual 15375 124.00
## Number of obs: 24, groups: UID, 6
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 665.50 60.73 4.00 10.959 0.000394
## poly(H2O2, 3)1 -2257.81 175.36 12.00 -12.876 2.2e-08
## poly(H2O2, 3)2 121.60 175.36 12.00 0.693 0.501220
## poly(H2O2, 3)3 -410.10 175.36 12.00 -2.339 0.037475
## genotypeYFP-ALC1 94.83 85.88 4.00 1.104 0.331430
## poly(H2O2, 3)1:genotypeYFP-ALC1 116.02 247.99 12.00 0.468 0.648277
## poly(H2O2, 3)2:genotypeYFP-ALC1 -477.25 247.99 12.00 -1.924 0.078329
## poly(H2O2, 3)3:genotypeYFP-ALC1 -170.66 247.99 12.00 -0.688 0.504421
##
## (Intercept) ***
## poly(H2O2, 3)1 ***
## poly(H2O2, 3)2
## poly(H2O2, 3)3 *
## genotypeYFP-ALC1
## poly(H2O2, 3)1:genotypeYFP-ALC1
## poly(H2O2, 3)2:genotypeYFP-ALC1 .
## poly(H2O2, 3)3:genotypeYFP-ALC1
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) pl(H2O2,3)1 pl(H2O2,3)2 pl(H2O2,3)3 gYFP-A p(H2O2,3)1:
## pl(H2O2,3)1 0.000
## pl(H2O2,3)2 0.000 0.000
## pl(H2O2,3)3 0.000 0.000 0.000
## gntYFP-ALC1 -0.707 0.000 0.000 0.000
## p(H2O2,3)1: 0.000 -0.707 0.000 0.000 0.000
## p(H2O2,3)2: 0.000 0.000 -0.707 0.000 0.000 0.000
## p(H2O2,3)3: 0.000 0.000 0.000 -0.707 0.000 0.000
## p(H2O2,3)2:
## pl(H2O2,3)1
## pl(H2O2,3)2
## pl(H2O2,3)3
## gntYFP-ALC1
## p(H2O2,3)1:
## p(H2O2,3)2:
## p(H2O2,3)3: 0.000
cat("AIC: ", AIC(fit12))
## AIC: 224.6933
simres <- simulateResiduals(fittedModel = fit12)
plot(simres)

Final Result
fit <- fit11
output <- coef(summary(fit))
output <- output[grep("H2O2", rownames(output)),]
rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype", paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1], sep = " in " )
# suggested result table
kable(output, row.names = T)
| H2O21 in WT |
-1.6341883 |
0.0513987 |
-31.7943497 |
0.0000000 |
| H2O22 in WT |
0.0884280 |
0.0513987 |
1.7204317 |
0.1046320 |
| H2O23 in WT |
-0.3093886 |
0.0513987 |
-6.0193863 |
0.0000179 |
| H2O21: WT vs. YFP-ALC1 |
0.0461736 |
0.0726887 |
0.6352232 |
0.5342627 |
| H2O22: WT vs. YFP-ALC1 |
-0.3507551 |
0.0726887 |
-4.8254386 |
0.0001864 |
| H2O23: WT vs. YFP-ALC1 |
-0.1157002 |
0.0726887 |
-1.5917212 |
0.1310090 |
write.table(output, file = "FigureS5N_Stats.txt", quote = F, sep = "\t", row.names = T, col.names = NA)